{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Load the required packages\n",
    "using ODE\n",
    "using JLD\n",
    "\n",
    "# Set precision to quadruple precision (QP)\n",
    "# QP has 113 bits\n",
    "set_bigfloat_precision(113);\n",
    "\n",
    "# Define the right-hand function\n",
    "function arenstorf(t,x)\n",
    "    \n",
    "    mu = parse(BigFloat,\"0.012277471\"); ms = 1 - mu;\n",
    "    r1 = vecnorm(x[1:2]-[-mu,0]);\n",
    "    r2 = vecnorm(x[1:2]-[ ms,0]);\n",
    "    \n",
    "    dx = Vector{BigFloat}(4);\n",
    "    dx[1] = x[3];\n",
    "    dx[2] = x[4];\n",
    "    dx[3] = x[1] + 2*x[4] - ms*(x[1]+mu)/r1^3 - mu*(x[1]-ms)/r2^3;\n",
    "    dx[4] = x[2] - 2*x[3] - ms*     x[2]/r1^3 - mu*     x[2]/r2^3;\n",
    "    return dx;\n",
    "end\n",
    "\n",
    "# Set up the initial conditions\n",
    "t0 = parse(BigFloat,\"0.0\");\n",
    "T = parse(BigFloat,\"17.0652165601579625588917206249\");\n",
    "x0=[parse(BigFloat,\"0.994\"),\n",
    "    parse(BigFloat,\"0.0\"),parse(BigFloat,\"0.0\"),\n",
    "    parse(BigFloat,\"-2.00158510637908252240537862224\")];\n",
    "\n",
    "# Set the tolerance as machine epsilon \n",
    "Tol = parse(BigFloat,\"1e-33\");\n",
    "\n",
    "# Solve and get the solution at T = tEnd\n",
    "(t,x_ref) = ode78(arenstorf,x0,[t0,T];\n",
    "reltol=Tol,abstol=Tol,points=:specified);\n",
    "\n",
    "# Save the solution to a file\n",
    "save(\"refSolAren.jld\",\"x_ref\",x_ref[2]);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.4.2",
   "language": "julia",
   "name": "julia-0.4"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.4.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
